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Abstract 


Spacecraft are being designed based on LS-DYNA simulations of water landing impacts. The Elemental 
Water Impact Test (EWIT) series was undertaken to assess the accuracy of LS-DYNA water impact 
simulations. Phase 1 of the EWIT series featured water impact tests of a 20-inch hemisphere dropped 
from heights of 5 feet and 10 feet. The hemisphere was outfitted with an accelerometer and three 
pressure gages. The focus of this report is the correlation of analytical models against test data. Three 
analytical models were used: 


1. Closed-Form Solution 

2. 2-D Axisymmetric LS-DYNA Model 

3. 3-D Quarter Symmetry LS-DYNA Model 


The closed form solution was found to over-predict the peak acceleration. The peak pressure from the 
closed-form solution is infinite; however, the solution does provide a reasonable prediction for the late - 
time pressure decay. 

The 2-D axisymmetric LS-DYNA model was evaluated for a broad range of modeling parameters. The 
two parameters that were found to be the most critical were mesh density and fluid- structure coupling 
stiffness. The 2-D model provided reasonable predictions for the acceleration histories. Pressure 
histories were poorly predicted; however, reasonable predictions were obtained for the impulse as 
determined by integrating the pressure history. Increasing mesh density resulted in better predictions for 
the late-time pressure decay, but did not improve the predictions of the early -time peak pressure. 

The 3-D quarter symmetry LS-DYNA model was used for further evaluations of mesh density and fluid 
structure coupling stiffness. The 3-D model produced predictions for the acceleration histories that 
better matched the test data than the 2-D model. The pressure histories from the 3-D model were equally 
as poor as the 2-D model. Again, the impulse was predicted far more accurately than the pressure. 

The following guidelines were proposed for defining the fluid mesh density and fluid- structure coupling 
stiffness for LS-DYNA simulations. These guidelines might not be feasible for spacecraft water landing 
design simulations due to limits in available computer resources. 

1. The element size of the fluid mesh should be no larger than l/200th of the radius of curvature of the 
structure. 

2. The coupling stiffness curve should be established based on the maximum expected coupling pressure 
at a penetration equal to 1/1 Oth the size of the fluid elements. 

Guidelines were also offered for judging the adequacy of the fluid- structure coupling definition. 

1. The coupling stiffness is too soft if the peak coupling pressure is toward the middle of the contact 
patch rather than at the perimeter. 

2. The coupling stiffness is too soft if the acceleration time history exhibits oscillations that can be 
related to the structure bouncing on the coupling surface. 

3. The coupling stiffness is too stiff if the coupling pressure distribution appears as a checkerboard 
pattern of isolated pressure spikes. 

4. The coupling stiffness is too stiff if the coupling pressure histories appear as a series of isolated 
spikes. 
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1. Introduction 


Spacecraft are being designed based on LS-DYNA [1] simulations of water landing impacts. The 
Elemental Water Impact Test (EWIT) series was undertaken to assess the accuracy of LS-DYNA water 
impact simulations. Phase 1 of the EWIT series featured water impact tests of a 20-inch hemisphere 
dropped from heights of 5 feet and 10 feet. 


2. Simulations 
2.1. LS-DYNA 

LS-DYNA is a general purpose transient dynamic finite element code capable of simulating complex real 
world problems. LS-DYNA’s strength is in the modeling of impact problems. An explicit time 
integration scheme is used in which there is no equilibrium check and iteration of the solution between 
time steps. This approach works only because the time step is restricted to be smaller than the shortest 
stress wave transit time for any element in the model. 

The modeling of contact between bodies in LS-DYNA is accomplished via a penalty method. Contact is 
detected when the nodes of one body pass through the face or edges of the elements of another body. 
Preloaded penalty springs are then inserted to push the bodies apart. One consequence of this approach is 
one body must always penetrate another body before contact is detected. Another consequence is that 
there is a finite contact stiffness at the interface between the bodies that is entirely nonphysical. 

LS-DYNA has a limited capability to model a fluid using Arbitrary Lagrangian-Eulerian (ALE) meshes. 
In the ALE approach, each time step begins with a mesh that is conceptually similar to the Lagrangian 
meshes used to model structures. LS-DYNA determines the deformation of the fluid that occurs during 
the time step, then moves, or advects, the mesh back to its original configuration and treats the fluid as 
having moved through the mesh. The result is that the nodes of the mesh do not move. Instead, the 
volume fraction of the fluid in each element is changed. The fluid in the ALE mesh can flow, compress, 
and impart momentum; however, the ALE mesh does not offer a full Navier-Stokes fluid flow solution. 

LS-DYNA features many parameters that can be adjusted in modeling of fluid-structure interaction 
problems. An overarching parameter is the mesh density. As with all finite element codes, a smaller 
element size generally produces results that are more accurate. In addition to mesh density, there are 
many other solution parameters that can be selected. Several of these are discussed in the following 
sections. 


2.2. ALE Solution Parameters 

The *CONTROL_ALE card includes several parameters, most of which should be allowed to take default 
values. The number of solution cycles per mesh advection step, NADV, should ideally be specified as 1, 
though it can be set to a higher value if necessary to reduce solution time. Various options for the mesh 
advection method, METH, are available. The developer recommends METH=2, which is second order 
accurate. Other options can be tried if necessary to reduce solution time. Several parameters are 
available for mesh advection smoothing. The developer warns that the smoothing algorithm is not based 
in physics. It is recommend that smoothing be turned off, ALAC=-1. 
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2.3. Element Section Parameters 

For a 2-D mesh, *SECTION_ALE2D offers two axisymmetric element formulations. ELFORM=14 is 
area weighted and ELFORM=15 is volume weighted. The code developer strongly recommends 
ELFORM=14. For a 3-D mesh, *SECTION_SOLID offers two treatments of an air space above a fluid, 
ELFORM=l 1 allows the air to be treated as a fluid with its own equation of state. ELFORM=12 treats 
the air as a void space. The code developer recommends ELFORM=l 1. 


2.4. Material Parameters 

Fluids with an equation of state utilize *MAT_NULL. The most important material property for the 
impact problem is the material density, RO. The tensile pressure cutoff, PC, is typically set to either zero 
or a very small negative number. If the air is treated as a void space (ELFORM=12), it may be desirable 
to set PC equal to negative one atmosphere (-14.7 psi) in order to mimic the fact that cavitation does not 
occur until the pressure drops to zero absolute pressure. The LS-DYNA Keyword User’s Manual [2] 
implies that a non-zero value must be specified in order for the fluid to be able to cavitate. 


2.5. Equation of State Parameters 

The two most commonly used equations of state for water are *EOS_GRUNEISEN and 
*EOS_LINEAR_POLYNOMIAL. *EOS_GRUNEISEN is required for water at high pressure as would 
occur in the simulation of an underwater explosion. For the range of pressures seen in the water impact 
problem, the two options provide similar relationships between pressure and volumetric strain. 
*EOS_LINEAR_POLYNOMIAL is readily usable for both water and air. 


2.6. Fluid Structure Coupling Parameters 

The *CONSTRAINED_LAGRANGE_IN_SOLID card features many options. The default for the 
number of fluid-structure coupling points on the surface of each structural element, NQUAD, is a 2 x 2 
array of coupling points. If the structural elements are larger than the fluid elements, a larger number of 
NQUAD points should be specified. Perhaps the most important coupling parameter is the coupling 
stiffness, PFAC. If PFAC is a positive value, the coupling stiffness is based on the solution time step and 
the mass of the nodes on either end of the penalty coupling springs. The default is PFAC=0.1. A 
negative value for PFAC points to a user-specified coupling stiffness curve. The curve specifies the 
coupling pressure as a function of the distance the fluid penetrates past the coupling surface. The user can 
also specify the fluid element minimum volume fraction at which coupling is activated, FRCMIN. The 
default is FRCMIN=0.5. Problems with excessive penetration of fluid into the structure can be addressed 
by reducing FRCMIN; however, there is no firm physical basis for this. Simulations with low values of 
FRCMIN have been observed to show the fluid moving away from the structure prematurely. 


2.7. Boundary Conditions 

The typical approach for modeling a water or soil block in an impact problem is to restrain the sides and 
bottom of the mesh in the directions normal to each face. The problem with this approach is that stress 
waves are reflected when they reach the mesh boundary. Usually, dispersion of the stress waves results in 
negligible pressure amplitude by the time the stress wave arrives back in the area of the impact. An 
alternative approach is to designate a layer of elements along the sides and bottom of the mesh as ambient 
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pressure reservoir elements. This is done via the ambient element type option, AET=4, on the 
*SECTION cards. The ambient pressure reservoir elements allow material to flow in and out of the mesh 
while holding the pressure constant. 


2.8. Gravity and Atmospheric Preload 

It is customary in water and soil impact problems to initialize gravity at lg instantaneously at the start of a 
simulation. The result is that stress waves oscillate through the material throughout the solution duration. 
This typically has negligible effect on the solution as the amplitude of the pressure oscillations is similar 
to the hydrostatic head at the base of the mesh, which is much less than the pressure due to the impact. A 
more accurate method is to ramp gravity to lg prior to the impact. The result is a stable hydrostatic stress 
state in the material provided that the ramp time is much longer than the time required for a stress wave to 
propagate through the depth of the mesh. If there is no atmospheric overpressure, the pressure in the 
material at the surface will be zero, which could allow material to cavitate prematurely. If the air is 
treated as a void (ELFORM=12), premature cavitation can be prevented by setting the tensile cutoff 
pressure on the *MAT_NULL input, PC, equal to negative one atmosphere (-14.7 psi). If the air is 
modeled as a fluid (ELFORM=ll), the *INITIAL_HYDROSTATIC_ALE and 
*ALE_AMBIENT_HYDROSTATIC options can be used to preload the air and water meshes with 
atmospheric and hydrostatic pressure. For the *INITIAF_HYDROSTATIC_AFE and 
*ALE_AMBIENT_HYDROSTATIC options, *EOS_FINEAR_POFYNOMIAF must be used to define 
the equation of state for the air and water. 


3. Tests 

3.1. Test Configuration 

The drop tests were performed in a 15 -foot diameter above-ground swimming pool located in the Room 
123 high-bay area of Building 1293 A at NASA Fangley Research Center (LaRC). The depth of the water 
in the pool was approximately four feet. The test pool was located inside a 24-foot diameter above- 
ground swimming pool to catch any over splash. A foam pad was placed under the liner at the bottom of 
the inner pool to cushion bottom impacts. The test article was suspended above the test pool on a boom 
extended from a forklift. A line hanging from the test article was used to measure the drop height. The 
test set-up is illustrated in Figure 1 . 
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Figure 1. Test Set-Up 


3.2. Test Article 

The test article was a hemispherical aluminum shell with an outside diameter of 20 inches. The general 
configuration of the test article is illustrated in Figure 2. The hemisphere had a shell thickness of 0.19 
inches and was filled with bismuth ballast to an approximate depth of 2.5 inches at the apex. Eight 
phenolic ribs served as internal attachment points for instrumentation, which included an accelerometer at 
the center and three pressure transducers at a distance of approximately 2.25 inches from the apex. The 
locations of the pressure transducers are illustrated in Figure 3. An aluminum lid with a thickness of 0.25 
inches served to keep water away from the instrumentation. The hemisphere was lifted into position for 
the drop test via three lift points at the perimeter of the lid. The weight of the hemisphere with 
instrumentation was 48 lb. 



Figure 2. Hemisphere Configuration 
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3.3. Test Results 

Three tests were performed at a 5 -foot drop height and three were performed at a 10-foot drop height. For 
comparison with simulation data, the acceleration histories were filtered using a 180 Hz Butterworth filter 
to remove the vibratory structural response. Figures 4 and 5 show raw and filtered acceleration histories 
for the 5-foot drops. Figures 6 and 7 show the acceleration histories for the 10-foot drops. An artifact of 
the filtering is a small negative dip in acceleration prior to impact. 



Figure 4. Raw Acceleration Histories for 5-foot Drops 
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Figure 5. Filtered Acceleration Histories for 5 -foot Drops 



Figure 6. Raw Acceleration Histories for 10-foot Drops 
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Figure 7. Filtered Acceleration Histories for 10 -foot Drops 

Histories from the pressure transducers for the 5-foot and 10-foot drops are shown in Figures 8 and 9. 

The sampling rate for both pressure and acceleration data was 50,000 Hz. The TDAS PRO data 
acquisition system used for the tests is equipped with an analog anti-aliasing filter with a cutoff frequency 
of 4300 Hz. The anti-aliasing filter can be bypassed, but it was used for these tests. The sampling rate 
was relatively high for this type of test; however, the duration of the pressure peaks is so short that even 
this relatively high sampling rate provides just a few points to define the pressure peaks. As a 
consequence, some of the pressure histories exhibit truncated peaks. Even for the pressure histories that 
do not show truncated peaks, it is inevitable that the data points sampled during the test did not capture 
the absolute pressure peak. 

One pressure history from each test that does not show a truncated peak was chosen for comparison 
against simulation data. The chosen pressure histories are Gage 3 of Test 1, Gage 3 of Test 2, and Gage 
2 of Test 3 for the five-foot drops, and Gage 1 of Test 4, Gage 2 of Test 5, and Gage 3 of Test 6 for the 
ten-foot drops. The peaks of the pressure histories do not occur simultaneously due to small angles of 
pitch that exist upon impact of the hemisphere with the water. In the figures, the relative timings of the 
peaks are as they occurred during the test. For comparisons between tests, the pressure histories have 
been arbitrarily shifted so that several peaks from different tests align. 
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Figure 8. Pressure Histories for 5 -foot Drops 
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Figure 9. Pressure Histories for 10 -foot Drops 


4. Simulation Models 
4.1. Closed-Form Solution 

A closed-form solution for the acceleration and pressure histories of a sphere impacting water was 
developed by A.P. Cappelli and J.P.D. Wilkinson [3]. The equations are as follows. 
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1/2 


2(i+r) 2 




c = (2 /? t) 1 / 2 

8V2p(RV 0 t) 3 / 2 
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3M 


Equation 1 
Equation 2 

Equation 3 
Equation 4 


Where: A = Acceleration 

P = Pressure 
p = Water Mass Density 
V 0 = Initial Velocity 
M = Mass 

R = Radius of Curvature at Impact Point 
r = Polar Distance from Impact Point 
c = Maximum Radius of Wetted Shell Surface 
t — Time after Impact 
Y = Non-Dimensional Parameter 

Cappelli and Wilkinson note problems with Equation 2. The predicted pressure is singular at the time 
that the contact patch radius, c, reaches the polar distance, r, which results in the equation producing no 
real answer for the highest pressure at a given polar distance. Cappelli and Wilkinson recommend 
dropping the second term of the pressure equation as it is small everywhere except for a negative 
singularity at the perimeter of the contact patch. Cappelli and Wilkinson further note that the remaining 
term of the pressure equation predicts nonphysical negative pressure in the late time. 


4.2. 2-D LS-DYNA Simulation 

A model of a two-dimensional axisymmetric slice of the penetromter was used to evaluate the sensitivity 
of the acceleration and pressure to several LS-DYNA parameters. The model configuration is illustrated 
in Figure 10. 
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The hemisphere was treated as linear elastic. Aluminum material properties were specified for the shell 
and the cover. The internal ballast was bismuth. The material density of the bismuth ballast was adjusted 
to achieve the desired total weight of 48 lb. The water was modeled with an equation of state. The air 
was treated as a vacuum. Gravity was ramped to lg during the first 0.01 seconds to obtain a stable 
hydrostatic pressure state in the water prior to impact. The initial height and velocity of the hemisphere 
were adjusted so that impact occurred at the desired velocity at 0.01 seconds. The LS-DYNA cards that 
describe the coupling stiffness, material properties, and initial conditions are provided in Appendix A. 

Three mesh variants of the model were used with element edge lengths of 0.1, 0.05, and 0.025 inches in 
the area of the initial contact for both the hemisphere shell and the water. The meshes are illustrated in 
Figure 1 1 . 
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Figure 11. Mesh Variants for Two-Dimensional Axisymmetric Model 


The pressure transducers used for the tests had a diameter of approximately one-eighth inch. For the 
simulations, coupling pressure was recorded using the DBFSI option in LS-DYNA. The pressure was 
averaged for the elements that span 0.1 inches at the pressure transducer location. This was one element 
for the 0.1 -inch mesh, two elements for the 0.05-inch mesh, and four elements for the 0.025-inch mesh. 

The baseline coupling stiffness was designated Curve 10. Curve 10 featured 100 psi at a penetration of 
0.005 inches. Simulations were also performed with variants of Curve 10 with the stiffness scaled by 
factors of 0.1, 10, and 100. 


4.3. 3-D LS-DYNA Simulation 

A three-dimensional quarter model of the hemisphere was used for further evaluations of coupling 
stiffness and mesh density. The three-dimensional model was considered necessary because the 2-D and 
3-D sections of the LS-DYNA code were written by different developers. The theory is the same for 
both, but there is no certainty that the implementation behaves in exactly the same manner, so there is no 
certainty that the findings concerning one modeling approach are directly applicable to the other. The 
three-dimensional model is illustrated in Figure 12. 
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Figure 12. Three-Dimensional Quarter Model 

The hemisphere was treated as rigid. Both the air and the water were modeled with equations of state. 
Gravity was initiated instantaneously and options within LS-DYNA were employed to initiate the air at 
one atmosphere and the water at one atmosphere plus hydrostatic pressure. The initial height of the 
hemisphere above the water was set so that impact occurred immediately after the start of the simulation. 
The bottom and outermost layers of elements of the water and air meshes were defined as reservoir 
elements. The reservoir elements allow material to flow in or out of the mesh in order to maintain 
constant pressure at the boundary. The LS-DYNA cards that describe the coupling stiffness, material 
properties, and initial conditions are provided in Appendix B. 

Three water meshes were used. The baseline mesh featured an element size of 0.1 inches near the initial 
impact. The overall mesh radius was 15 inches, the air height was 1.2 inches, and the water depth was 
13.2 inches. The first variant featured an element size of 0.05 inches near the initial impact and the same 
overall mesh radius and depth as the baseline mesh. The second variant was created by halving the 
dimensions of the 0.05-inch mesh. This resulted in an element size of 0.025 inches near the initial impact 
and an overall mesh radius of 7.5 inches, an air height of 0.6 inches, and a water depth of 6.6 inches. The 
0.025-inch element model is illustrated in Figure 13. Due to the limited extent of the mesh, the 0.025- 
inch element model was suitable only for short duration simulations to determine early -time impact 
response. 
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Figure 13. 0.025 -inch Fluid Element Mesh Variant 

For all simulations, the element size for the hemisphere was 0.1 inches, which is similar to the size of the 
pressure transducer used in the tests. Pressures from the simulations were recorded for a single element 
using the DBF SI option. 

The baseline coupling stiffness was designated Curve 12. Curve 12 featured 100 psi at a penetration of 
0.01 inches. This was half the stiffness of the baseline Curve 10 used for the 2-D simulations. The 
rationale for halving the coupling stiffness for the 3-D simulations was that the baseline 3-D mesh 
featured elements with edge lengths twice as large as the baseline 2-D mesh. For the baseline mesh with 
an element size of 0.1 inches, simulations were performed with variants of Curve 12 that were scaled by 
factors of 0.1 and 10 for the stiffness. Simulations were also performed with the coupling stiffness 
default, PFAC = 0.1. LS-DYNA determines the default coupling stiffness based on the stiffness required 
to achieve a vibration period equal to the critical time step size for solution stability. 


5. Test and Simulation Correlation 
5.1. Closed-Form Results 

Accelerations for the closed-form solution for the 5 -foot and 10-foot drops are plotted with test data in 
Figures 14 and 15. The test data was filtered using a 180-Hz Butterworth filter to eliminate vibratory 
structural response not represented in the analytical models. The closed-form solution exhibits no noise 
and requires no filtering. The results show that the closed-form solution substantially over predicted the 
peak acceleration. This is partly due to the closed-form solution not taking into account the reduction in 
the velocity that occurs as a consequence of the deceleration during the early stages of the impact. For the 
5-foot drop, the velocity would have reduced from 215 in/sec at impact to 168 in/sec at 0.01 seconds after 
impact. For the 10-foot drop, the velocity would have reduced from 304 in/sec to 206 in/sec. 


13 



.jest l 

Test 2 

-Test 3 


Closed-Form 

Solution 


Figure 14. Closed-Form Solution Acceleration History for 5-foot Drop 
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Figure 15. Closed-Form Solution Acceleration History for 10-foot Drop 

The closed-form solution pressure histories are plotted with test data in Figures 16 and 17. The closed- 
form solution is asymptotic at the perimeter of the contact patch, but provides a reasonable prediction of 
the pressure decay. 
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Figure 16. Closed-Form Solution Pressure History for 5 -foot Drop 
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Figure 17. Closed-Form Solution Pressure History for 10-foot Drop 


5.2. 2-D LS-DYNA Results 

Acceleration histories for the three mesh variants are compared to test accelerometer histories in Figures 
18 through 20 for the 5-foot drops and Figures 21 through 23 for the 10-foot drops. All acceleration 
histories were filtered using a 1 80-Hz Butterworth filter and were time-shifted to facilitate comparison. 
The figures for the 5-foot drops show close agreement with the test results. For the 10-foot drops, the 
simulations over-predicted the acceleration peak. Both the 5 -foot and 10-foot drops show the simulations 
had a much slower acceleration decay rate in the late time, after approximately 0.008 seconds. 
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Figure 18. Acceleration Histories for 5-foot Drop for 0.025-inch Element Size 
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Figure 19. Acceleration Histories for 5-foot Drop for 0.05-inch Element Size 
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Figure 20. Acceleration Histories for 5 -foot Drop for 0.1 -inch Element Size 
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Figure 21. Acceleration Histories for 10-foot Drop for 0.025-inch Element Size 
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Figure 22. Acceleration Histories for 10-foot Drop for 0.5-inch Element Size 
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Figure 23. Acceleration Histories for 10-foot Drop for 0.1 -inch Element Size 

The acceleration histories show that the filtered peak acceleration is relatively insensitive to both the 
mesh density and the coupling stiffness for the range of element sizes and coupling stiffnesses considered. 
This is illustrated in the bar charts shown in Figure 24 for the 5-foot drops and Figure 25 for the 10-foot 
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drops. The tests produced filtered peaks of approximately 1 lg for the 5 -foot drops and 22g for the 10- 
foot drops. 
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Figure 24. Peak Accelerations for Mesh and Coupling Stiffness Variants for 5 -foot Drops 
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Figure 25. Peak Acceleration for Mesh and Coupling Stiffness Variants for 10-foot Drops 

Pressure histories were extracted from the LS-DYNA simulations using the DBFSI option, which 
translates the coupling forces into an interface pressure. The DBFSI pressure does not necessarily match 
the pressure seen in the fluid. This is illustrated in Figure 26 for a 5-foot drop of the 0.025-inch mesh 
model with the 10 x Curve 10 coupling stiffness. The DBFSI pressure reflects the actual loading on the 
coupling surface whereas the fluid pressure represents the pressure at the center of the elements of the 
fluid mesh. The relative mesh densities for the fluid and structure along with the density of the fluid 
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coupling points on the structure (NQUAD) all factor into how closely the pressure in the fluid matches 
the interface pressure. Further mesh refinement would lead to better agreement between the pressure in 
the fluid and the pressure at the coupling surface. 
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Figure 26. Pressure in Fluid versus DBFSI Pressure at Coupling Surface 

Pressure histories for the three mesh variants are compared to test pressure histories in Figures 27 through 
29 for the 5-foot drops and Figures 30 through 32 for the 10-foot drops. The pressure histories are 
unfiltered. Pressure histories for the 10 x Curve 10 simulations are shown in Figure 33 for the 5-foot 
drops and Figure 34 for the 10-foot drops. The simulation pressure histories have been averaged over an 
edge length of 0.1 inches, which approximates the pressure transducer size of 0.125 inches. The 
simulation results show that refining the mesh results in better agreement with the test data for the late 
time decay of the pressure history; however, there is no convergence in matching the peak pressure of the 
test data. 


20 



Test 1 Gage 3 

— Test 2 Gage 3 

— — — Test 3 Gage 2 

0.025-inch Mesh, Curve 10x0.1 

0.025-inch Mesh, Curve 10x1.0 

0.025-inch Mesh, Curve 10x10 

0.025-inch Mesh, Curve 10x100 


Figure 27. Pressure Histories for 5-foot Drop for 0.025-inch Element Size 
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Figure 28. Pressure Histories for 5-foot Drop for 0.05-inch Element Size 
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Figure 29. Pressure Histories for 5-foot Drop for 0.1-inch Element Size 
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Figure 30. Pressure Histories for 10-foot Drop for 0.025-inch Element Size 
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Figure 31. Pressure Histories for 10-foot Drop for 0.05-inch Element Size 
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Figure 32. Pressure Histories for 10-foot Drop for 0.1 -inch Element Size 
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Figure 33. Pressure Histories for 5 -foot Drop for 10 x Curve 10 Mesh Variants 
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Figure 34. Pressure Histories for 10 -foot Drop for 10 x Curve 10 Mesh Variants 

The pressure histories show that the peak pressure is highly sensitive to both the mesh density and the 
coupling stiffness. The results also show that this sensitivity is reduced as the mesh is refined; however, 
this does not mean that the peak pressures from the simulations converge toward the peak pressures from 
the test. This is illustrated in the bar charts shown in Figure 35 for the 5-foot drops and Figure 36 for the 
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10-foot drops. The tests produced peak pressures of approximately 120 psi for the 5 -foot drops and 200 
psi for the 10-foot drops. 
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Figure 36. Peak Pressure for Mesh and Coupling Stiffness Variants for 10 -foot Drops 

Figures 37 and 38 show the impulse at the pressure transducer locations from the LS-DYNA simulations. 
The impulse was calculated by integrating the pressure history for the first 0.01 seconds following impact. 
Despite the drastic differences in the shapes of the pressure histories, the impulse is relatively insensitive 
to the mesh density and the coupling stiffness. The impulse varies between 0.069 psi-sec and 0.073 psi- 
sec for the simulations of the 5-foot drops and between 0.112 psi-sec and 0.117 psi-sec for the simulations 
of the 10-foot drops. This explains why the acceleration histories were very similar despite drastic 
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differences in the pressure histories. For simulation of the response of the overall structure, the exact 
profile of the pressure pulse is probably not important as long as the impulse is correct; however, the 
profile of the pressure pulse may be important for the design of small penetrations or crushable materials. 

A similar computation was performed for the pressure histories from the tests. The impulse for the 5 -foot 
drops was 0.116 psi-sec for Gage 3 of Test 1 and 0.0.084 psi-sec for Gage 3 of Test 2. The impulse for 
the 10-foot drops was 0.122 psi-sec for Gage 2 of Test 5 and 0.119 psi-sec for Gage 3 of Test 6. The 
pressure histories for Gage 2 of Test 3 and Gage 1 of Test 4 were not long enough to determine the 
impulse for the full period. 


u 

m 


01 

tfli 



■ 0.1-inch Mesh 

■ 0.05-inch Mesh 

■ 0.025-inch Mesh 


Figure 37. Impulse for Mesh and Coupling Stiffness Variants for 5 -foot Drops 
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Figure 38. Impulse for Mesh and Coupling Stiffness Variants for 10 -foot Drops 
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In addition to the effect the coupling stiffness has on the peak of the pressure history, it also affects the 
shape of the pressure contours within the fluid. Also, a coupling stiffness that is too soft can result in 
substantial penetration of fluid through the fluid-structure interface. This is illustrated in Figure 39 for the 
simulation of a 5-foot drop of the 0.05-inch mesh. The 0.1 x Curve 10 coupling stiffness produces 
smooth pressure contours but shows significant penetration of the fluid into the structure. The 100 x 
Curve 10 coupling stiffness shows no penetration, but the pressure contours show multiple isolated 
patches of high pressure. For the 100 x Curve 10 coupling stiffness, the impact occurs as a series of 
localized hits rather than as one continuous impact. 
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Figure 39. Pressure Contours for Coupling Stiffness Variants for 5-foot Drop 


Sensitivity studies were performed for several simulation parameters. The parameters varied are listed in 
Table 1. The baseline simulation for the study was the 0.05-inch element model dropped from 5 feet. An 
exception was made for the divot model, which was based on the 0.025-inch element model in order to 
provide a finer mesh for defining the divot. Bar charts illustrating the peak filtered acceleration and peak 
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pressure are shown in Figures 40 and 41 . Most of these parameters had relatively little effect on the peak 
filtered acceleration and peak pressure. 


Table 1. Parameters Varied for Sensitivity Study 


Parameter 

Description 

Baseline 

ELFORM 

Element formulation for 2-D 
axisymmetric mesh. 

14 

EOSGRUNEISEN 

Form of the equation of state for 
the water. 

EOSLINEARPOLYNOMIAL 

FRCMIN 

Volume fraction at which a fluid 
element is considered to be 
coupled to the structure. 

0.5 

METH 

Method for advecting fluid mesh 
between solution time steps. 

3 

PC 

Tensile cutoff pressure for fluid. 

14.7 psi 

Instantaneous Gravity 

Method of applying gravity. 

Gravity ramped over 0.01 seconds 
prior to impact. 

Rigid Model 

Treating the hemisphere model 
as rigid rather than flexible. 

Flexible model. 

61.3 lb Model 

Weight of hemisphere. 

48.0 lb 

Divot 

Adding a divot with a depth of 
0.02 inches in way of the 
pressure transducer location. 

No divot. 
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Figure 40. Peak Filtered Accelerations for Sensitivity Study Variants 
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Figure 41. Peak Pressures for Sensitivity Study Variants 


Regarding these parameters, the following should be noted. 


ELFORM=15 does not compute the stress correctly and should not be used according to 
the LS-DYNA developers at LSTC. 

The value of FRCMIN should be changed with caution as a low value can result in the 
fluid mesh moving before the impacting object arrives. 

METH=2 is described in the FS-DYNA Keyword User’s Manual [2] as “second order 
accurate”, whereas METH=1 and METH=3 (Baseline) are described as “first order 
accurate”. 

PC=-14.7 is used to mimic the overpressure of the atmosphere. The result is that the 
water at the surface of the fluid mesh does not cavitate until a tensile pressure of 14.7 psi 
is reached. 

Initiating gravity instantaneously results in an oscillating pressure field in the fluid 
whereas ramping gravity results in a stable hydrostatic pressure field. 

The measured weight of the hemisphere was 48 lbs. The 61.3 lb weight was based on the 
basic geometry assuming that the bismuth was sold rather than machined away to 
accommodate instrumentation. 

Increasing the weight of the hemisphere would be expected to reduce the peak deceleration. This is 
consistent with the simulation results. The reduced deceleration would result in a modestly higher 
velocity at the time that the pressure pulse reaches the transducer location, which would be expected to 
result in a higher pressure peak. The simulation showed negligible change in the peak pressure. The 
effect of the increased weight on the pressure history was so small that no clear effect was perceptible in 
the simulation. 

The idea of adding a divot to the simulation model stemmed from the difficulty of obtaining a perfectly 
flush mounting for a pressure transducer in the test article. For the 2-D axisymmetric model, the 
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representation of a recessed transducer mounting as a divot was not exactly correct as the divot actually 
represented an axisymmetric groove; however, the effect was expected to be similar. The 0.02 -inch divot 
resulted in a significant spike in the pressure as the water passed over it. This is illustrated in Figure 42 
for a 5-foot drop of the 0.025-inch mesh with the 10 x Curve 10 coupling stiffness. This was consistent 
with findings from experiments with pressure transducers that showed that slightly recessing the pressure 
transducer resulted in a higher pressure reading. 
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Figure 42. Pressure Pulse Caused by a 0.02-inch Divot 


The conclusion from the results from the 2-D model was that the most important parameters for the study 
of test versus analysis correlation are mesh density and coupling stiffness. Other parameters either have 
little effect or lack any physical basis for making a change. As a consequence, mesh density and coupling 
stiffness were the focus of further studies performed with the three-dimensional model of the hemisphere. 


5.3. 3-D LS-DYNA Results 

A series of coupling stiffness variants for the baseline mesh (0.1 -inch elements) were analyzed. In 
Figures 43 and 44, the acceleration histories for the 5 -foot and 10-foot drops are compared to test data. 
For both drop heights, the simulation acceleration histories correlate well with the test data. The variant 
with the lowest coupling stiffness (0.1 x Curve 12) shows oscillatory behavior, which resulted from the 
hemisphere bouncing on the coupling surface. 
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Figure 43. Acceleration History for 5 -foot Drop Simulations 
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Figure 44. Acceleration History for 10-foot Drop Simulations 
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The pressure histories from the simulations are compared to test data in Figures 45 and 46. The 
simulation pressure histories show a strong sensitivity to coupling stiffness. The baseline coupling 
stiffness (Curve 12x1) shows the closest match to the peak pressures measured in the test; however, the 
shape of the curve from the simulation does not match the shape of the curves from the test. 
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Figure 45. Pressure Histories for 5-foot Drop Simulations 
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Figure 46. Pressure Histories for 10-foot Drop Simulations 

For simulations of the 5-foot drops, the impulse during the first 0.01 seconds after impact varied between 
0.052 psi-sec for the Curve 12 x 10 variant and 0.067 psi-sec for the Curve 12x1 variant. For the 5-foot 
drop test data, the range was 0.084 psi-sec to 0.116 psi-sec. For simulations of the 10-foot drops, the 
impulse varied between .098 for the Curve 12 x 10 variant and 0.100 for the PFAC Default variant. For 
the 10-foot drop test data, the range was 0.119 psi-sec to 0.122 psi-sec. 

Figures 47 and 48 illustrate the distribution of the coupling pressure in the simulations of the 5-foot and 
10-foot drops. The softest coupling stiffness (0.1 x Curve 12) produced a broad pressure patch without a 
well-defined “Coliseum Effect” ring of high pressure at its perimeter. The highest coupling stiffness 
(PFAC Default) produced a scattering of high pressure spikes rather than continuous rings of pressure. 
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Figure 47. Pressure Distribution at 0.001 seconds after Impact for 5-foot Drop Simulations 
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Figure 48. Pressure Distribution at 0.001 seconds after Impact for 10 -foot Drop Simulations 

Simulations were run for fluid mesh variants with element sizes of 0.05 inches and 0.025 inches. The 
element size of the hemisphere was 0.1 inches for all simulation variants. The 0.05-inch fluid mesh was 
created by subdividing the 0.1 -inch elements of the baseline fluid mesh. The 0.025-inch fluid mesh was 
created by halving the dimensions of the 0.05-inch mesh. As a consequence, the 0.025-inch fluid mesh 
has just half the overall extent of the 0.1 -inch and 0.05-inch fluid meshes. For the 0.025-inch fluid mesh, 
the hemisphere began to interact with the mesh boundary at approximately 0.002 seconds. This limited 
the 0.025-inch fluid mesh simulation to the very early time response. 

For the 0.05-inch fluid mesh, the coupling stiffness was made twice as stiff as the Curve 12 baseline. For 
the 0.025-inch fluid mesh, coupling stiffness was made four times as stiff as the Curve 12 baseline. The 
increase in coupling stiffness for the finer meshes was motivated by the belief that the maximum 
allowable penetration distance should be a function of the fluid element size rather than an absolute 
dimension. Also, the number of fluid coupling points (NQUAD) was increased from 2 for the 0.1 -inch 
baseline fluid mesh to 3 for the 0.05-inch fluid mesh and 6 for the 0.025-inch fluid mesh. This was 
necessary to avoid the smaller fluid elements slipping between the coupling points of the 0.1 -inch 
elements of the hemisphere. 
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Acceleration histories for the 5 -foot and 10-foot drops for the fluid mesh variants are compared to test 
data in Figures 49 and 50. The 0.025-inch fluid mesh is not included in the acceleration comparison 
because the hemisphere begins interacting with the fluid mesh boundary before the peak acceleration was 
reached. The results show that the refinement from a 0.1 -inch mesh to a 0.05-inch mesh makes negligible 
difference in the acceleration histories. 



Test 1 

Test 2 

-Test3 

0.1" Mesh Curve 12 x 1 
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Figure 49. Acceleration Histories for Mesh Variants for 5 -foot Simulations 
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Figure 50. Acceleration Histories for Mesh Variants for 10-foot Simulations 

Pressure histories for the fluid mesh variants are compared to test data in Figures 51 and 52. The pressure 
histories show that the baseline 0.1 -inch mesh offers the best prediction for the peak pressure, but does 
not follow the pressure decay curves from the tests and has a second peak that is not physically 
reasonable. Refinement to a 0.05-inch mesh reduces the pressure peak below the value from the tests and 
results in a pressure history that oscillates about the pressure decay curves from the tests. Further 
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refinement to a 0.025-inch mesh further reduces the peak pressure and results in a pressure history that 
closely follows the pressure decay curves from the tests. 
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Figure 51. Pressure Histories for Mesh Variants for 5-foot Simulations 
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Figure 52. Pressure Histories for Mesh Variants for 10 -foot Simulations 

The simulation results for the 3-D model had peak filtered accelerations that were very similar to the test 
results. For the 5-foot drops, the simulations under-predicted the peak accelerations by less than 15%. 

For the 10-foot drops, the simulations modestly over-predicted the peak accelerations. The acceleration 
results support the use of a 15% allowance for design for loads that are proportional to the acceleration. 
The same cannot be concluded for the pressure histories. The peak pressures can be either over-predicted 
or under-predicted based on the mesh density and the coupling stiffness. 
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6. Establishing Simulation Parameters without Test Data 
6.1. Width of Pressure Pulse 

A series of tests were performed with additional pressure transducers mounted in the hemisphere. The 
time required for the pressure pulse to propagate from one transducer location to another provides an 
indication of the speed at which the pressure pulse travels. This can be used in conjunction with the 
observed time duration of the pressure pulse to estimate of the width of the pressure pulse. Figure 53 
shows the arrangement of pressure transducers. Table 2 lists the coordinates of the pressure transducers. 



Figure 53. Pressure Transducer Arrangement 


Table 2. Pressure Transducer Locations 


Sensor ID 

X (in) 

Y (in) 

MSI 

0.8125 

2.125 

MS2 

1.250 

-3.1875 

MS3 

0.875 

-2.0625 

PCB1 

-0.875 

-2.0625 

PCB2 

3.0625 

0.875 

PCB3 

-1.750 

-4.1875 


The histories for the pressure transducers for a 10-foot drop are shown in Figure 54. The peak pressure 
decreases, and the time duration increases, as the pressure pulse moves away from the apex. Table 3 
provides an estimate of the width of the pressure pulse based on the time duration of the pulse and the 
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velocity of the pulse as determined from the arrival time at subsequent transducer locations. The 
calculation suggests a pulse width of approximately 0.2 inches, which is approximately l/50 th of the 
radius of the hemisphere. 


Drop #2: Pressure Transducers 



Figure 54. Pressure Histories for Multiple Pressure Transducers 
Table 3. Pressure Pulse Width Estimate 


Gage 

Arrival 

Time 

(sec) 

X (in) 

Y (in) 

Radius 

(in) 

Velocity 

(in/sec) 

Pulse 

Width 

(sec) 

Pulse 

Width 

(in) 

PCB1 

4.711 

-0.875 

-2.0625 

2.2404 

986.17 

0.0002 

0.197 

MS2 

4.7122 

-1.25 

-3.1875 

3.4238 

857.41 

0.0002 

0.171 

PCB3 

4.7135 

-1.75 

-4.1875 

4.5385 

- 

- 

- 


6.2. Establishing Mesh Density and Coupling Stiffness 

The ability of the model to adequately simulate the pressure pulse is dependent on having a mesh density 
that adequately describes the shape of the pressure pulse. If it is assumed that a minimum of four 
elements are needed across the width of the pressure pulse, and it is estimated that the width of the 
pressure pulse is approximately l/50 th of the radius of curvature of the structure, then the fluid element 
size should be no larger than l/200 th of the radius of curvature. This leads to a fluid element size of 0.05 
inches for the hemisphere. This is the element size for the middle mesh variant for the 2-D axisymmetric 
model. It is half the element size of the baseline 3-D quarter model. 
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The approach recommended by LSTC in the LS-DYNA Keyword User’s Manual [2] for defining the 
coupling stiffness curve is as follows. 

“The curve consists of {0,0} as the first point and {maximum allowable penetration (MAP), estimated 
maximum coupling pressure (EMCP)} as a second point. MAP may be a small penetration with respect to 
the minimum ALE element width (maybe 10% or less). EMCP can be estimated from a maximum fluid 
pressure observed from a previous run when leakage first occurs. This curve may be scaled to vary the 
stiffness of the coupling spring. The approach is to gradually increase the coupling stiffness until leakage 
stops. The best coupling stiffness is one which provides just enough force to prevent leakage and not 
more ” 

LSTC subsequently amended this recommendation in e-mail and verbal communications. The amended 
recommendation is to base EMCP on the maximum pressure observed in a simulation with the PFAC 
default (PFAC = 0.1). There remains ambiguity as to whether this is the peak pressure seen anywhere in 
the model or at a presumed pressure transducer location and whether it should be the pressure seen in the 
fluid or at the fluid-structure interface. FS-DYNA does permit the specification of different fluid- 
structure coupling definitions (*CONSTRAINED_FAGRANGE_IN_SOFID) for different regions of the 
model; however, the coupling definition cannot be changed during the simulation. 

For the hemisphere, the recommended mesh density and coupling stiffness are described in Table 4. 


Table 4. Recommendations for Element Size and Coupling Stiffness 


Parameter 

Criteria 

Value 

Radius of Curvature, R 

Known Dimension 

10 in 

Fluid Element Size, F 

L = R/200 

0.05 in 

Maximum Allowable 
Penetration, MAP 

MAP = L/10 

0.005 in 

Estimated Maximum Coupling 
Pressure, EMCP 

From LS-DYNA simulation 
with PFAC = 0.1 

100 psi 


The maximum pressure is expected to be more than 10 psi but less than 1000 psi, so 100 psi was chosen 
as an order of magnitude estimate for EMCP. This resulted in a curve that is identical to Curve 10 (100 
psi at 0.005 inch) used for the 2-D axisymmetric model. It is twice as stiff as Curve 12 (100 psi at 0.01 
inch) used for the 3-D quarter model. It is recognized that even with these guidelines, the peak pressure is 
not likely to be predicted with a high degree of accuracy. 


6.3. Evaluating Model Performance without Test Data 

When test data is unavailable, the analyst must make judgments on the performance of the fluid-structure 
interaction to determine whether the coupling stiffness is too high or too low. Through the course of these 
and other simulations, several aspects of the model performance have been observed that provide 
indications regarding whether the behavior of a simulation is physically reasonable. 

One way to assess the adequacy of the fluid-structure interaction is to observe the degree of penetration 
and leakage through the structure. The following distinction is made between penetration and leakage. 
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Penetration is fluid that passes through the coupling interface but is still resisted by the 
coupling force. 

Leakage is fluid that escapes through the coupling interface and is not resisted by the 
coupling force. 

Figure 55 illustrates the distinction between penetration and leakage. 


O.l-inch Structure Mesh 
with 5x5 Array of 
Fluid CouplingPoints 


Outer area of contact 
patch shows 
penetration of fluid 
through the coupling 
surface; however, the 
fluid is still resisted by 
the couplingstiffness. 



0.025-inch Fluid Mesh 


Inner area of contact 
patch shows leakage of 
fluid that has escaped 
through the coupling 
surface and is no longer 
resisted by the coupling 
stiffness. 


Figure 55. Description of Penetration and Leakage 


Since LS-DYNA utilizes a penalty method in which springs are added to resist the passage of material 
through the coupling interface, some degree of penetration must always occur in order for there to be a 
coupling force. Penetration can be reduced by increasing the coupling stiffness (PFAC). 

Leakage should not occur at a properly functioning fluid-structure interface. Leakage may occur if the 
relative mesh size is such that not every fluid element is resisted by coupling points on the structural 
elements. Leakage can be alleviated by increasing the density of the coupling points (NQUAD). 

Two symptoms of an excessively soft coupling stiffness were observed. The first was a coupling pressure 
pattern that lacks a well-defined “Coliseum Effect” ring of high pressure at the outer perimeter of the 
contact patch. The second was low frequency oscillation of the acceleration history as a consequence of 
the oscillatory system created by the mass of the structure and the coupling stiffness. These symptoms 
are illustrated in Figure 56. 


41 



Early time pressure 
distribution lacks a 
well-defined 
"Coliseum Effect" 
ring at the 
perimeter of the 
contact patch. 


Acceleration history 
exhibits oscillations 
due to ringing of the 
contact force. 



Test 1 

- - - Test 2 

Test 3 

Curve 12x0.1 

Curve 12x1 

Curve 12 xlO 

PFAC Default 


Figure 56. Symptoms of an Excessively Soft Coupling Stiffness 

Two symptoms were also observed for an excessively stiff coupling stiffness. The first was an interface 
coupling pressure pattern that appears as a set of isolated point loads. The second was a pressure history 
that appears as a series of isolated high pressure spikes. These symptoms are illustrated in Figure 57. 
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Pressure 

distribution appears 
as a series of 
isolated point loads. 


Pressure history 
appears as a series 
of isolated spikes. 
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0.1" Mesh Curve 12x0.1 
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0.1" Mesh Curve 12 x 10 

0.1” Mesh PFAC Default 


Figure 57. Symptoms of Excessively High Coupling Stiffness 

These criteria can be used to evaluate coupling stiffnesses that are at the extremes. Unfortunately, this 
still leaves a broad range of coupling stiffnesses that are either too stiff or too soft but give results that are 
not obviously wrong. 


7. Conclusions and Recommendations 

The results for the hemisphere model demonstrate that LS-DYNA can accurately predict acceleration 
histories for a broad range of mesh densities and coupling stiffnesses. The prediction of pressure histories 
is a more difficult problem. LS-DYNA predictions of pressure histories improve as the fluid mesh is 
refined; however, the required mesh size for a real world problem may be impractical due to large model 
size requiring excessive disk space and unacceptably long run times. 
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The simulation results for the 3-D hemisphere model had peak filtered accelerations that were very 
similar to the test results. For the 10-foot drops, the simulations modestly over-predicted the peak 
accelerations. For the 5-foot drops, the simulations under-predicted the peak accelerations by less than 
15%. This supports the use of a 15% allowance for design for loads that are proportional to the 
acceleration. The same cannot be concluded for the pressure histories. The peak pressures can be either 
over-predicted or under predicted based on the mesh density and the coupling stiffness. 

In the absence of test data, the analyst must make a judgment regarding whether the fluid-structure 
coupling interface is performing adequately. The following guidelines are proposed. 


1. The element size of the fluid mesh should be no larger than l/200 th of the radius of 
curvature of the structure. 

2. The coupling stiffness curve should be established based on the maximum expected 
coupling pressure at a penetration equal to 1/1 0 th the size of the fluid elements. 

It is recognized that the above guidelines may result in an impractical model size with impractical mn 
time and that the pressure histories are not likely to be realistic even if the guidelines are followed. As a 
consequence, the following guidelines forjudging the adequacy of the performance of the coupling 
surface are offered. 


1. The coupling stiffness is too soft if the peak coupling pressure is toward the middle of the 
contact patch rather than at the perimeter. 

2. The coupling stiffness is too soft if the acceleration time history exhibits oscillations that 
can be related to the structure bouncing on the coupling surface. 

3. The coupling stiffness is too stiff if the coupling pressure distribution appears as a 
checkerboard pattern of isolated pressure spikes. 

4. The coupling stiffness is too stiff if the coupling pressure histories appear as a series of 
isolated spikes. 
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Appendix A: Two-Dimensional Axisymmetric Model 


The following are the LS-DYNA cards that control the material properties, contact, and initial conditions 
for the axisymmetric model. These particular cards are for the 0.05-inch mesh for the 10-foot drop with 
the Curve 10 x 1 coupling stiffness. Node numbers and the number of set segments for pressure output 
differ depending on the mesh density. 


* KEYWORD 800000000 
$ IMPACT AT_0 . 01 SECONDS 
$ RAMPED GRAVITY 
$ ELFORM=14 FOR AIR AND WATER 
$ VARIABLE D3PLOT TIME STEP 
''TITLE 

2-D Axisymmetric Model of 20" Hemisphere 
*CONTROL_TERMI NATION 

$ ENDTIM ENDCYC DTMIN ENDENG ENDMAS 

0.0200 0.0 0.0 0.0 
$*CONTROL_PARALLEL 
$4 

*CONTROL_HOURGLASS 

1 0.1 
*CONTROL_ENERGY 


$# HGEN RWEN 

SLNTEN 

RYLEN 





2 2 







*CONTROL_OUTPUT 
$# NPOPT NEECHO 

NREFUP 

IACCOP 

OPIFS 

IPNINT 

IKEDIT 

IFLUSH 


1 3 
$# IPRTF 

0 

$ 

$*DATABASE_RBDOUT 
$ 0.0000020 
*DATABASE_NODOUT 
0.0000020 

*DATABASE_HISTORY_NODE 

$ idl id2 i d 3 id4 i d 5 id6 id7 id8 

4000020 4000024 4001553 4001557 4001975 

* DATAB AS E_B I NARY_D 3 P LOT 
$# dt lcdt beam npltc 

0.0000200 7778 

*DEFINE_CURVE 

7778 0 1.0 1.0 0.0 0.0 0 

0.00000 0.002000 

0.00990 0.002000 

0.01000 0.000020 

0.01210 0.000020 

0.01212 0.002000 

0.02000 0.002000 

* DATAB AS E_B I NARY_D3THDT 
999.999 

$ 

* S ET_M U LTI - MAT E R I A L_G ROU P_L I ST 
$ 3=Vacuum 2=Soil 
123 

2 

*ALE_MULTI-MATERIAL_GROUP 
2 1 

3 1 

$ 

* S ET_PART_L I ST 
1 

2 3 

*SET_PART_LIST 
11 

1 11 12 21 22 

*CONTROL_ALE 


$ 

DCT 

3 

NADV 

1 

METH 

3 

AFAC 

-1.0 

BFAC 

CFAC 

DFAC 

EFAC 

$ 

START 

END 

AAFAC 

VFACT 

PRIT 

EBC 

PREF 

NSIDBEC 
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$ 


’'CONSTRAINED 

_L AG RAN G E_I N_SO L I D_ 

.EDGE 





$# slave 

master 

sstyp 

mstyp 

nquad 

ctype 

di rec 

mcoup 

11 


1 

0 

0 

3 

4 

2 

-123 

$# start 


end 

pfac 

f ri c 

f rcmi n 

norm 

normtyp 

damp 

0 


0 

-10 

0 

0 

0 

0 

0 

$# cq 


hmi n 

hmax 

i 1 eak 

pi eak 

lcidpor 

nvent 

i block 

0 


0 

0 

0 

0.1 




$# i box id 

i penchk 

i ntforc 

ialesof 

1 agmul 

pfacmm 

thkf 


0 

c 


0 

1 

0 

0 




-p 

*DEFINE_CURVE 








$ lcid 


si dr 

sfa 

sfo 





10 



1.0 

1.00 







0.000 


0.0 







0.005 


100.0 





-p 

*SECTION_SHELL 








$ SID 


ELFORM 

SHRF 

NIP 

PROPT 

IRID 

ICOMP 

SETYP 

1 


14 

0.83333 





1 

$ Tl 


T2 

T3 

T4 

NLOC 

MAREA 

IDOF 

EDGSET 

1 . 


1 . 

1 . 

1 . 





*SECTION_SHELL 








$ SID 


ELFORM 

SHRF 

NIP 

PROPT 

IRID 

ICOMP 

SETYP 

98 


14 

0.83333 





1 

$ Tl 


T2 

T3 

T4 

NLOC 

MAREA 

IDOF 

EDGSET 

1 . 


1 . 

1 . 

1 . 





*SECTION_SHELL 








$ SID 


ELFORM 

SHRF 

NIP 

PROPT 

IRID 

ICOMP 

SETYP 

99 


14 

0.83333 





1 

$ Tl 


T2 

T3 

T4 

NLOC 

MAREA 

IDOF 

EDGSET 

1 . 


1 . 

1 . 

1 . 





*SECTION_ALE2D 








$ SID 

ALEFORM 

AET 

ELFORM 





2 


11 


14 





*SECTION_ALE2D 








$ SID 

ALEFORM 

AET 

ELFORM 





3 

c 


11 


14 





4> 

* PART 









vacuum 









$ PID 


SECID 

MID 

EOSID 

HGID 

GRAV 

ADAPT 

TMID 

2 


3 

3 

0 

3 

0 



4> 

* PART 









water 









$ PID 


SECID 

MID 

EOSID 

HGID 

GRAV 

ADAPT 

TMID 

3 

c 


2 

2 

2 

2 

0 



4> 

* PART 









Hemi sphere 









$ PID 


SECID 

MID 

EOSID 

HGID 

GRAV 

ADAPT 

TMID 

1 


1 

1 

0 

1 

0 



SET SEGMENT 

11 

ELEMENT 







11 


1 

1 

0 

1 

0 



SET SEGMENT 

12 

ELEMENT 







12 


1 

1 

0 

1 

0 



SET SEGMENT 

21 

ELEMENT 







21 


1 

1 

0 

1 

0 



SET SEGMENT 

22 

ELEMENT 







22 


1 

1 

0 

1 

0 



Bi smuth 









98 


98 

98 

0 

1 

0 



Cover 









99 

$ 


99 

99 

0 

1 

0 



$ 

*MAT_VACUUM 









$ mid 


rho 







3 1 

.116E-11 







*MAT_NULL 









$ mid 


ro 

PC 

mu 

terod 

cerod 

ym 

pr 

2 9 

. 3365e-5 

-14.7 

1.6300E-7 0. 

.0000000 0, 

.0000000 



*MAT_ELASTIC 








1 

O.i 

000253 

10 . 2e6 

0.33 



0 
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4.6E6 


0.33 


0 


*MAT_ELASTIC 

98 0.0006279 

*MAT_ELASTIC 

99 0.000253 10.2E6 0.33 


$ 

*EOS_LINEAR_POLYNOMIAL 

2 0.0000000 3 . 11574e5 0.0000000 


0.0 

1.0 



$ 




’’'Hourglass 




$ HGID 

IHQ 

QM 


3 

1 

l.E-6 


2 

1 

l.E-6 


1 

1 

0.1 


* S ET_NOD E_L I ST_G EN E RATE 



$ Sid 




111 




$ blbeg 

blend 



4000001 

4100000 



* INITIAL_VELOCITY 



$ nsid 




111 




$ vx 

vy 

vz 

vxr 

0 . 

-302.48 

0 . 


’’'LOAD BODY Y 




1 

386.1 



’'DEFINE CURVE 




1 

0 

1.0 

1.0 


0.0 


0.0 


0.01 


1.0 


100.0 


1.0 


$ 

’’'INCLUDE 

2d_hemi sphere_water_block_rev4_0p050. k 
*DE FIN E_TRANS FORMATION 
100 

$ al a2 a3 

TRANSL 0.0 13.04 0. 

$ 

*INCLUDE_TRANSFORM 

2d_20i nch_hemi sphere_0p050. k 


$ 

i dnoff 

ideoff 

idpoff 

i dmof f 


4000000 

4000000 


$ 

idroff 




$ 

fctmas 

fcttim 

fctlen 

fcttem 

$ 

t ran id 





100 





$ 

*END 


0 

0.0000000 0.0000000 0.0000000 


vyr 

vzr 


0.0 

0.0 

0 


a4 

a5 

a6 

a7 

off 

i dfoff 

iddoff 



i ncoutl 
1 


47 



Appendix B: Three-Dimensional Quarter Model 


The following are the LS-DYNA cards that control the material properties, contact, and initial conditions 
for the three-dimensional quarter model. These particular cards are for the 0.1-inch mesh for the 10- 
foot drop with the Curve 12 x 1 coupling stiffness. For the 0.05-inch mesh, the number of fluid-structure 
coupling points, NQUAD, is increased to 3. For the 0.025-inch mesh, NQUAD is increased to 6. 


’'KEYWORD 800000000 
*TITLE 

3D QUARTER MODEL OF HEMISPHERE 
*CONTROL_TERMI NATION 

$ ENDTIM ENDCYC DTMIN ENDENG ENDMAS 

0.02 0.0 0.0 0.0 

*CONTROL_HOURGLASS 

1 0.1 

*CONTROL_ENERGY 

$# HGEN RWEN SLNTEN RYLEN 

2 2 

*CONTROL_OUTPUT 

S# NPOPT NEECHO NREFUP IACCOP OPIFS IPNINT IKEDIT IFLUSH 

1 3 

$# IPRTF 

0 

$ 

* DATAB AS E_G L STAT 
0.0000020 
* DATAB AS E_MATSUM 
0.0000020 
*DATABASE_RBDOUT 
0.0000020 

* DATAB AS E_B I NARY_D 3 P LOT 
$# dt lcdt beam npltc 

0.0002 

*DATABASE_BINARY_FSIFOR 
$ dt 

$ 0.0000200 7778 

0.0002 

* DATAB AS E_B I N ARY_D 3TH DT 
999.999 

$ 

*SET_PART_LIST 

1 

2 3 22 23 

* S ET_PART_L I ST 

222 

2 22 

* S ET_PART_L I ST 

323 

3 23 

*ALE_MULTI-MATERIAL_GROUP 
$ sid idtype 

222 0 

323 0 

* S ET_M U LTI - MAT E RAI L_G ROU P_L I ST 

123 

2 

*CONTROI ALE 


$# 

dct 

nadv 

meth 

afac 

bfac 

cfac 

dfac 

efac 


2 

1 

2 

- 1.0 





$# 

start 

end 

aafac 

vfact 

prit 

ebc 

pref 

nsidebc 

$ 







14.7 


^CONSTRAINED. 

.LAGRANGE. 

_IN_SOLID 






$# 

slave 

master 

sstyp 

mstyp 

nquad 

ctype 

di rec 

mcoup 


1 

1 

1 

0 

2 

4 

2 

-123 

$# 

start 

end 

pfac 

f ri c 

f rcmi n 

norm 

normtyp 

damp 


0 

0 

-12 

0 

0 

0 

0 

0 

$# 

cq 

hmi n 

hmax 

ileak 

pleak 

lcidpor 

nvent 

i block 


0 

0 

0 

0 

0.1 




$# 

i box id 

i penchk 

i ntforc 

ialesof 

1 agmul 

pfacmm 

thkf 
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0 


0 

1 

0 

0 



$ 








*DEFINE_CURVE 







$ lend 


si dr 

sfa 

sfo 




12 



1.0 

1.00 






0.000 


0.0 






0.010 


100.0 




$ 








*SECTION_SHELL 






$ sid 


el form 

shrf 

ni p 

propt 

i rid 

i comp 

1 


1 0. 

,8333333 





$ tl 


t2 

t3 

t4 

nloc 

marea 

idof 

0.1875 


0.1875 

0.1875 

0.1875 




*SECTION_SOLID 






$ SID 


ELFORM 

AET 





2 


11 






*SECTION_SOLID 






$ SID 


ELFORM 

AET 





5 


11 






*SECTION_SOLID 






$ SID 


ELFORM 

AET 





22 


11 

4 





*SECTION_SOLID 






$ SID 


ELFORM 

AET 





25 


11 

4 





$ 








* PART 








Hemi sphere 
$ PID 


SECID 

MID 

EOSID 

HGID 

GRAV 

ADAPT 

1 


1 

1 





$ 








* PART 








vacuum 
$ PID 


SECID 

MID 

EOSID 

HGID 

GRAV 

ADAPT 

2 


2 

2 

2 

2 

0 


$ 








* PART 








water 

$ PID 


SECID 

MID 

EOSID 

HGID 

GRAV 

ADAPT 

3 


3 

3 

3 

3 

0 


$ 








* PART 








vacuum 
$ PID 


SECID 

MID 

EOSID 

HGID 

GRAV 

ADAPT 

22 


22 

2 

2 

2 

0 


$ 








* PART 








water 

$ PID 


SECID 

MID 

EOSID 

HGID 

GRAV 

ADAPT 

23 


23 

3 

3 

3 

0 


$ 

*MAT_NULL 
$ mid 


rho 

PC 

mu 

terod 

cerod 

ym 

2 

1. 

127E-7 

0.0 




*MAT_NULL 
$ mid 


rho 

PC 

mu 

terod 

cerod 

ym 

3 

9 . 3365e-5 

0.0 1. 

, 6300E-7 0. 

.0000000 0. 

.0000000 

*MAT_RIGID 
$ mid 


ro 

e 

pr 




1 

0. 

001075 

le+07 

0.3 




$ emo 


coni 

con2 





1 


5 

7 





$ al 


a2 

a3 

vl 

v2 

v3 


$ 








$*EOS_IDEAL. 

_GAS 






$$ eosid 


cv 

cp 

cl 

c2 

to 

vO 

$ 2 


6.179E5 

8.651E5 

0.0 

0.0 

527.67 

1.0 

*EOS_LINEAR_POLYNOMIAL 






$ eosid 


cO 

cl 

c2 

c3 

c4 

c5 

2 


0.0 

0.0 

0.0 

0.0 

0.4 

0.4 

$ eO 


vO 






36.74 


0.0 







$ 

*EOS_LINEAR_POLYNOMIAL 


setyp 

edgset 


TMID 

TMID 

TMID 

TMID 

TMID 

pr 

pr 


c6 

0.0 
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cO cl c2 c3 c4 c5 c6 

14.7 3 . 11574e5 0.0000000 0.0000000 0.0000000 0.0000000 

vO 

0.0 


IHQ 

1 


$ eosid 

3 

$ eO 

0.0 

$ 

''Hourglass 
$ HGID 

2 

3 1 

* S ET_NOD E_L I ST_G E N E RATE 
$ Sid 

111 

$ blbeg blend 

4000001 4900000 

*INITIAL_VELOCITY 


$ nsid 


111 


$ vx 

vy 

304.41 

0. 

*LOAD_BODY_X 


1 

-386.1 

*DEFINE_CURVE 

1 

0 


0.0 

100.0 

*SET_PART_LIST 

$ sid 


781 


$ pi dl 

pid2 

2 

3 

* I N ITI AL_H YDROSTATIC_AL E 

$ SID 

SIDTYPE 

781 

0 

$ NID 

MMGBELOW 

103067 

1 

44 

2 

* S ET_PART_L I ST 

$ sid 


782 


$ pi dl 

pid2 

22 

23 

*ALE_AMBIENT. 

.HYDROSTATIC 

$ SID 

SIDTYPE 

782 

0 

$ NID 

MMGBELOW 

103067 

1 

44 

2 


*DEFINE_VECTOR 
$ Vi d Xt 

789 0. 

$ 

’'INCLUDE 

water_ai r_quarter_Opl. k 
*DE FIN E.TRANS FORMATION 
100 

$ al 

TRANSL -10.01 

$ 

*INCLUDE_TRANSFORM 
hemi sphere_mesh_revO . k 
$ idnoff ideoff 

4000000 4000000 

$ idroff 


QM 

l.E-6 

l.E-6 


vz vxr vyr 

0 . 


1.0 1.0 0.0 
1.0 
1.0 


VECID GRAVITY PBASE 

789 386.1 14.7 


VECID GRAVITY PBASE 

789 386.1 14.7 


yt zt xh 

0 . 0 . 1 . 


a2 a3 a4 

0.0 0.0 


idpoff idmoff idsoff 


vzr 


0.0 0 


yh zh cid 

0 . 0 . 


a5 a6 a 7 


idfoff iddoff 


$ fctmas fcttim fctlen fcttem incoutl 

1 

$ tranid 

100 

$ 

*END 
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